“Calhoun 


Institutional Archive of the Naval Postgraduate School 





Calhoun: The NPS Institutional Archive 
DSpace Repository 


Theses and Dissertations 1. Thesis and Dissertation Collection, all items 


1971 


A computer study of the buckling of an 
imperfect cylindrical shell 


Dantone, Joseph John 


Monterey, California ; Naval Postgraduate School 
http://ndl.handle.net/10945/15625 


This publication is a work of the U.S. Government as defined in Title 17, United 
States Code, Section 101. Copyright protection is not available for this work in the 
United States. 


Downloaded from NPS Archive: Calhoun 


Calhoun is the Naval Postgraduate School's public access digital repository for 
| (8 D U DLEY research materials and institutional publications created by the NPS community. 
«ist sia Calhoun is named for Professor of Mathematics Guy K. Calhoun, NPS's first 


NY KNOX appointed — and published -- scholarly author. 

ia) LIBRARY Dudley Knox Library / Naval Postgraduate School 

411 Dyer Road / 1 University Circle 
Monterey, California USA 93943 





http://www.nps.edu/library 


A COMPUTER STUDY OF THE 
BUCKLING OF AN IMPERFECT 
CYLINDRICAL SHELL 


KSZOLELSESED 


JOSEPH JOHN DANTONE 














aa 








a 





United States 
Naval Postgraduate Schoo! 






RY 
NAN AVAL POSTGRADUATE SCHOOL 
UV ONTEREY, CALIF. 93940 


= oa 










ANGOMEURER STUDY OF THE BUCKLING OF AN 
IMPERFECT CYLINDRICAL SHELL 











Joseph John Dantone 


Hiesis Aavisor: 


September 1971 


= oa 


Approved for public aclease; distribution unlimcted. 


ATCoOnpucer Study Of Ehebuckling of an 
Imperfect Cylindrical Shell 


by 


Joseph John Dantone 
Lieutenant, United! States Navy 
B.S., United States Naval Academy, 1964 


Submitted in partial fulfiilment of the 
requirements for the degree of 


MASTER OF SCIENCE IN AERONAUTICAL ENGINEERING 


from the 


NAVAL POSTGRADUATE SCHOOL 
September 1971 





ABSTRACT 


An existing computer program for the geometrically non- 
linear analysis of arbitrarily loaded shells of revolution 
was modified to incorporate initial geometric imperfections 
imeo thesbuckling solution of an axially loaded circular 
cylindrical shell. Several different boundary conditions 
were considered, and the predicted buckling loads were 
examined to determine the significance of the initial imper- 
fections and the boundary conditions for the shell stability 
problem. The computed buckling loads were compared to those 
found experimentally by Arbocz and Babcock at the California 
Pattee Of Technelegy;,. The competed cosy 
generally higher than the experimental results. It was 
indicated that imperfection sensitivity is dependent on the 


boundary conditions of the cylinder. 
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Tro. ENERODUCTION 


tee study of the instability of thin, isetropic, circular 
cylindrical shells has long been troubled with a lack of 
correlation between the theoretically predicted buckling 
loads and the experimental results; the experimental values 
of buckling are generally much lower than the theoretical 
solutions. This lack of agreement has been attributed to 
imperfection sensitivity, i.e., inherent imoerfections 
in the shape of the shell which are unaccounted for in the 
classical buckling theory trigger the buckling modes pre- 
maturely, and to the use of incorrect boundary conditions 
in the analysis, i.@., tife boundary conditions used in the 
analysis are not equivalent to those of the test specimen [1]. 

The recent development of several sophisticated digital 
computer programs for the geometrically nonlinear analysis of 
shell structures has made it possible to include both the 
initial imperfections and arbitrary boundary conditions in 
the analysis. For example, a digital computer program for 
the geometrical nonlinear analysis of arbitrarily loaded 
shells of revolution has been developed by Ball [2]. A 
modification has recently been made to this program to 
include the initial imperfections in the analysis. However, 
this modification does not eliminate all of the previous 
analytical difficulties because in order to predict the 


buckling load of a particular shell, the imperfect shape of 





that shell must be known, and the actual boundary conditions 
of the specimen must be known. With regard to the first 
problem, Arbocz and Babcock [3] have developed a means of 
mapping the topography of cylindrical shells and can thus 
describe the deviations from a perfect circular cylinder. 
The objective of this thesis is to use the imperfection data 
for one of the cylinders presented in Ref. [3] and Ball's 
modified computer program to obtain stability solutions for 
several boundary conditions. The buckling loads are compared 
with the experimental results presented in Ref. [3] and are 
compared with each other to determine the effect of the 


various boundary conditions. 
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It. DESCRIP PON er tHe COMPUTER PROGRAM 


The computer program used in this study is capable of 
analyzing thin isotropic shells under the following condi- 
Ehens: 

1) The geometric and material properties must be axi- 
symmetric, but may vary along a meridian. 

2) The initial imperfections, pressures, and tempera- 
tures must be symmetric about a meridional plane. 

3) The boundaries may be free, fixed, or elastically 
restrained. 

The governing partial differential equations are based 
Gn Senders” nOnlinear Shelli theory for che conditienief 
small strains and moderately small rotations (Appendix A). 
Inplane and normal forces are accounted for but rotary 
inertial forces are neglected. The set of partial differ- 
ential equations is reduced to an infinite number of sets 
of four second-order differential equations in the meridional 
and time coordinates by expansion of the dependent variables 
in a sine or cosine series in the argument n9 where @ is 
the circumferential coordinate and n is the Fourier index. 
Nonlinear coupling terms are treated as pseudo loads and 
trigonometric identities are used to uncouple the sets of 
equations. Derivatives with respect to the meridional 
coordinates are replaced by the conventional finite differ- 


ence approximation. This leads to sets of algebraic 





(n) (n) (n) 


equations in the four dependent variables U pa , W : 


(n) 


S 


and M where U, V, and W are the displacements in the 
meridional, circumferential, and normal directions respec- 
tively, and M. is the meridional bending moment per unit 
length. At each load or time step, an estimate of the 
solution is obtained by extropolation from solutions at 
the previous steps. The sets of equations are repeatedly 
solved using Gaussian elimination, and the pseudo loads 
are recomputed until the solution converges. If convergence 
is not achieved in a specified number of iterations, the 
load is reduced hy a factor of five. If, after a specified 
number of load reductions, the solution fails to converge 
the seeEren is terminated. If the load-displacement 
Beewior of tire shell Gs of the softening ¢ype, DHueckiims 
1s presumed to have occurred at the final load. 
The original program described in [2] was not designed 
to solve shell stability problems where initial imperfections 
had to be accounted for. However, it soon became apparent 
that such a capability was very desirable and the program 
was modified accordingly. The equations and FORTRAN state- 


ments that include the imperfection terms are presented in 


Appendix A. 
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III. BOUNDARY CONDITIONS 


Pieeme analysis, the boundary conditions for the cylin- 
der are on Ne or U, N og OF ave 8 or W, and dW/ds or Me 
where N. and N56 are membrane forces per unit length, and 
On 1s the transverse force per unit length. The boundary 


conditions are represented by the matrix equation 


Na U 
eG V 
[2] o tt [A] = 
S W 
oW 
os ae 


where [2] and [Aj are 4 x 4 matrices and {2} is al x 4 
column matrix. There is one such equation for each boundary 
of the shell. Note that the analysis assumes that the 
boundary conditions are the same for all values of n, where 
nis the mode number in the circumferential direction. 

The testing configuration used in Ref. [3] for the 
axially loaded shell is shown in Fig. 1. The ends of the 
shell are attached to an end ring and a load cell with 
Cerrelow, a low melting point alloy. The ends of the shell 


are clamped such that 
W = dW/ds = v= 0 ) 


The axial load on the test specimen was applied by means 


of four screws located at the corners of the square end 


ie 





Cerrolow 


Load Coll 


Devcon 


End Se 





FIG.1 CYLINDRICAL SHELL TESTING CONFIGURATION 


ie 





plates. At each load increment, the screws were adjusted 
to give as uniform a circumferential load distribution as 
possible. 

In order to investigate the probiem completely, several 
boundary conditions were formulated for the computer program. 
These formulations follow below. 


If the axial load is axisymmetric, then on each boundary 


wf9) 2 y a el re 
S O S 

end (2) 
ui) x 9g 7" SO oe 


This condition will be referred to as boundary condition A. 
Moce in Eq. (2) that the condition on Ne 1S a £LUNCELON Of ane 
This necessitates changes in the program to account for 
different values of the first element of 2 for n = 0 and 

n 70. These changes are given in Appendix B. 


Another possibility is 


wy 69) = N ; N = 0 ee oes 
S O S 


at one end of the cylinder, and 


0, Ne 0 nie OO. ee 
at the other end. This condition will be referred to as 
boumigery condition B. 
A third interpretation of the test set-up consists of 
rigid end plates incrementally shortening the cylinder. 


This is represented at the boundaries by 





(n) 4h ee 


c 
o 
T 
c: 
c 
| 
© 
= 
| 


at one end, and 


u'™ 2 9 = Oily eee 


at the other end. This condition will be referred to as 


beum@ary condition C. 
For purposes of comparison, a simply supported cylinder 
will also be analyzed using boundary conditions A, B, and 


C. For this case, equation (1) becomes 


W= M = V= 0 (3) 
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IV. UNSEAT Tee REECTIONS 


Arbocz and Babcock's topography mapping procedure is 
described in detail in Ref. [3]. The final imperfection 
data consisted of deflections from a perfectly circular 
cylinder, denoted here by W, around fifteen circumferences 
spaced 0.5 in. apart. The first and the last circumferen- 
tial mappings occurred at a distance of 0.125 in. from the 
ends of the shell. Forty-nine data points were taken around 
each circumference with the first and the last being the 
same point on one arbitrarily chosen meridian. 

As noted previously, Ball's program requires a meridional 
plane about wnich the imperfect cylinder is symmetric. This 
requirement is not generally met by imperfect cylinders; 
thus symmetry must be created in some optimum way. An 
investigation to determine the degree of symmetry was con- 
ducted by generating a cosine Fourier series for each 
circumference using every meridian of data as a base 
meridian. In all cases, the phase angles for each circum- 
ference varied considerably, thus eliminating the possibility 
that the imperfections were symmetric. Consequently, a 
least squares method was used to determine the best meridian 
of data points about which to create symmetry. This best 


meridian was used as the origin to form the series 


7 Ze ey 
W= 2 W cos (n@) 
n=0 
using the data points from ®@ = (Q to 0 = T. 





An attempt was made to create a finer meridional mesh 
of the imperfection data by passing a general polynomial 
through the fifteen known coefficients for each value of n 
in order to interpolate imperfection coefficients at inter- 
mediate stations. However, a polynomial of sufficient 
degree to pass through all of the stations resulted in 
unreliable coefficients at intermediate stations. Con- 
versely, a least squares fit using a lower order polynomial 
resulted in a reasonable coefficient representation, but 
the coefficients at the original fifteen stations were 
compromised to an extent that depended on the degree of 
the polynomial. Similar problems were encountered using 
a Fourier series representation. As a result of the dif- 
ficuity in establishing coefficients at intermediate stations, 
the fifteen original locations were selected as the stations 
for the numerical shell analysis and the conventional central 
finite difference scheme was used to obtain the slope of 
the imperfection shape away from the boundaries. At the 
boundaries the simplest forward and backward differencing 


schemes were used. 
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V. AN EIGENVALUE SOLUTION 


To provide a check on the validity of the solution 
from the computer program, buckling of a simply supported 
perfectly circular cylindrical shell was formulated as an 
eigenvalue problem and solved analytically. The set of 


uncoupled field equations can be given in the form 


me), cp)" h 4 pel] ce)" b 4 pe ™] 27 = tel 4 Ca) 


where 


ip 
7 = fy) yin wim) mi} 


on) 


The eanations for the values of the elements in the [E 
cr f™))], ig'™], and fe) matrices are given in Appendix C. 
Ber n= 0, the linear membrane solution v0 S pe Is 
assumed. Thus the buckling modes of the shell for the 
boundary conditions given in equations (2) (boundary condi- 


tion A) and (3) are 


(n) 


U = Cy cos mis/L 
y (P) = C. sin mvs/L 

(5) 
wf?) = C4 Sin mis/L 


0) a 
M. = Cy Sin mis/L 


where n # O, Cy - Cy are arbitrary constants, m denotes the 


axial mode number and L is the length of the cylinder. 
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Substituting equations (5) into eguations (4) results in an 


eigenvalue problem of the form 
[A] {tc} = -N  [B]{c} 
where [A] and [B]) are 4 x 4 matrices defined in Appendix C 
and 
{c} = {c,,C,,C,,C,}° 
ae oe 4 
Because B is Singular, the inverse formulation 


> [Al{c} = -(B){c} (6) 


ine 





Vi. SEELL DESCRIPTION AND INPUT DATA 


The shell used in this study was cylinder A-8 of Ref. 3 


Its dimensions and properties are: 


Radius, r,Ra = 4.0 in. 

minha mores 5 Itt 
iertenGloeen: = 0.00464 in. 

a 5 yaya” 115//cte. 
Poisson's = Bae 


ratio,v 


Six harmenics of the initial imperfection data were used 
as input to the program. They are n = 0, 9, 10, 11, 12, and 
13. The selection of these particular harmonics was based 
on the eigenvalue solution for the simply supported shell 
that predicated a minimum buckling load when m = 1 and n= 9, 
and the experimental observation in Ref. [3] that n = 13 
appeared to be the critical mode where m and n are the mode 
shapes in the meridional and circumferential directions 
respectively. The number of meridional stations used was 
fifteen because no reliable interpolation of the initial 
imperfections could be developed, as discussed above. The 
length of the shell, L,was taken as 8 in. with the first and 
the last data point stations assumed to be the stations at 
the edges of the shell. Convergence criterion was taken 


as 0.01. 


We 





As presented in Appendix A, the imperfection input data 


Comersts of the rotations . and a5” where 
k k 
= 0 (n) 7 ao 
k O k k O 


in which - 1s a reference stress level, and Ey 1S a 
reference elastic modulus, a is a reference length, and k 
denotes the station number. The data are introduced into 
the program through the subroutine IMPERF. The FORTRAN 
statements for IMPERF are given in Appendix D. The imper- 
fection data was provided to the author by the writers of 
no! inethe form of a deck of cards With the data W/h. 
In order to obtain a buckling load that was a very 
close approximation to the bifurcation icad or the elgen- 
value problem, the imperfections were reduced by a factor 
#0? This technique has been used successfully with this 


program by Stillwell [4] who used very small asymmetries 


in the load as the triggering mechanism. 


20 





VII. RESULTS AND DISCUSSION 


The solution to the eigenvalue problem given by equa- 
tion (6) resulted in a minimum buckling load for the simply 
supported cylinder of 53.02 lb/in. when m= 1 andn = 9. 

The solution from the computer program for the minute imper- 
Feetiens resulted in a buckling load of 54.0 lb/in for 
boundary condition A, simply supported edges, and displace- 
ment in harmonics n = 0 and n= 9. Some difference between 
the two results is to be expected because the analytical 
solution neglects the pre-buckling deformation due to edge 
constraint, and the numerical solution used only fifteen 
stations. 

The computer results for the simply supported cylinder 
with boundary conditions A, B, and C, with both minute and 
full imperfections,* are presented in Table I. These results 
show a relationship between the imperfection sensitivity of 
the shell and the boundary conditions. For example, boundary 
Gondition A shows a 15% reduction in the load carrying capa- 
city when the full imperfection is introduced. Boundary 
Gomm@rtion B shows a 10% reduction, while condition C reduces 
Ene Gapacity by 42.7%. 

The computer results for the clamped shell are given in 
Table II. The experimental buckling load of shell A-8 was 


* = 
menute = full amperfections x 10 o. 
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Tabi 2 


oe iy soUPLOREE Dp. CYLINDER 


a ee 
————— 


Loading Degree of Buckling Load 
Condition Imperfection (lb/in) 

A minute 54.0 

A afi Sid ig | 45.9 

B minute 69.0 

B eek 62:49 

C minute 7 O16 

& fee 45.0 

TABLE II 


Cramen lp CYLENDER 











Loading Degree of Buckling Load 
Condition Imperfection (ils / Bia) 

A minute i 

no 

a Sud convergence 

B minute 43.87 

B Uo 46.9 

@: minute B57 

G gS lll 49.1 
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homme tO be 32.87 lb/in [3]... From Table II: 

1) The buckling load for condition A with minute imper- 
fections is 7.5 lb/in.* 

Ze NO converged solution weould be found. for condition 
A with full imperfections. 

3) The buckling load for condition B is slightly higher 
with full imperfections than for the minute case. 

4) The buckling load for loading condition C is 
decreased by 42.8% when the full imperfections are 
matroduced. 

Thus, the stability of the clamped shell under condi- 
tion Ais drastically reduced below that of condition B and 
C. This is very puzzling and much effort was expended in 
Seeren cfi,en explanation 


uUnitontyinately:- none weg found. 
The computed result for condition C, although generally 
higher than the empirical result, shows the expected load 
reduction due to the initial imperfections [3]. Some of 
the discrepancy in magnitude could be due to the use of an 
effective length that was 0.25 in.less than the actual 
length of the shell since the greater length would lower 
the computed buckling load. The remainder may be due to 
the fact that the actual boundary condition is some combin- 
ation of conditions A and C. Another consideration is that 
the symmetric representation of the actual imperfections 


may also degrade the validity of the solution. 


* 
For this boundary condition, the imperfections were 
Geduced by a factor of 1079. 
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Tae circumferential harmenictm = 10 was found to exhibit 


(0) 0) 


the most rapid growth. Figure 2 displays Ne 


vs W 


boundary condition C. The sudden decrease in slope occur- 


(0) 


Bingeaaist Prior to reaching re. 


illustrates the buckling 


(non-convergence) phenomenon. 
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VII “CORCLUSIONS 


The imperfection sensitivity of a cylinder is dependent 
on the boundary conditions. 

The computer program yields more realistic results for 
the simply supported boundary condition than for the 
ellamped conditi@ime 

The initial imperfections result in a significant 
reqwermon in load Carrying Capacity in the simply 


supported case but less reduction than expected [3]. 


Further investigation is necessary to determine: 


(ie cae tie comerees ct 
experimental results? 
Why is the clamped boundary condition A so sensitive to 
imperfections in the computer program? 

Why does loading condition B, clamped, appear to he 
strengthened by the introduction of full imperfections 


in the computer solution? 





APPENDIX A 


This section presents the appropriate equations for 
including initial imperfections in the computer program. 


The strain-displacement relationships used in this analysis 


are 
=e = Ute ee 4 Ope 
S S S 
f. Wee : 2 2 
Eg = V UG pe U/R + W/R + (6°49 eee (Al) 
= ' ° = ' 
San = CVU Jr i ey fe a 95) 2 
where 
o =-0' + Us 
S a Se 
%, =-W/r + V/Ro (A2) 
co) = (V' + r'V/r - U’/r)/2 


and where Eor Ege represent Fourier coefficients of the 


E39 


reference surface strains, $9, Pa. and ae represent refer- 


6 
ence surface rotations, and Ro and R. are the principal 
radii of curvature. The superscript primes and dots denote 
partial differentiation with respect to s and @ respectively. 
Consider an imperfect cylinder. Before loading, the 
imperfections appear as a displacement from the perfect 
cylinder. Note that this is a displacement with no accom- 


panying strain. Call this displacement W, and the total 


displacement from the perfect shape W*. Then 


2g 





W* = W+ W 


where W is that part of the total displacement caused by 
the applied load and is associated with strain. Replacing 


W with W* in equations (Al) and (A2) and requiring that 


« @ os ae a 0 when U = V = W = O, leads to 
e. = U'+W/R_ + [ (-W'4U/R_)* - 2W' (-W'+U/R_) + 9°] /2 
Cg = Vi /rtr'U/r +W/Rp + [(-WYr + V/R oe 2" /x (-W* /x+V/R,) 
t o71/2 | 
een U7 re V/r + (-WEsU7R.) (-W /raV7R W' 
ewer +P.) - W Ae Gwe 7R_))/2 
Besititing tif imperfection zetations 4_ and 9, as 
. = —W' and oy = -W'/r 


equations (A2) become 


20 a= 2 
= iy +W/R. + Cn 2a 2 ly 2 


= 
Ss 
Eee = "Ve7e +r u/r + WR, + (07 +25 wy) + 6%) /2 
G 6 6 6° 6 
a, +U°/r - r'V/r + Ae Pah D200) 72 
From Ref. [2] 
Oo. \2]/ @ oO co 
jee = »} 6!” cos ne|* = = g!) cos née 
2 Oo n=0 o n=0 








Redefine B such that 
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Oo co 
= 1} a cos n@ = 0° + 20.0. 
o n=0 
me | , (n) ~ (n) (n) 
= 5 LZ @es no); +2! £ 9 cos né Lae cos né 
S S 
Oo n=0 n=0 n=0 
Simaelarly 
a © 5 
E z By cos n@ = dot 20,0, 
o n=0 
C co 9 co ( co 
— oy ocos nee 2} MY) sin ng 2 ee Sin né@ 
E © = 8 6 
fe) n=1 | geal = 
and 
a co = _ 
z 2 B gSin n@ = PE 6 + ae + 9%, 
Oo n=1 








n=l 
where 

=(n) _ _ dwi?? 

w= ~~ de 
and 

gm) z sony, 

G 0 
= ——(n) 

pice = si/a, 9 = r/a, and wy (P) = ee 7 (ac0) 


The changes to the program are in SUBROUTINE PHIBET 


given below. 
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wee 


~~ 


L)*PHT( J) +PHIXC JS) *PHT CT) 
*PHT(JS)+PHIX(J)*PHT(T) 


—_ 


hwy OO ee” 


MAXS(1G) sMAXSY(10) ,18010,10) ,J8(16,10),10(10 
+GA=V2(M) +FEN®U2Z(0IM) *RRA)¥25 


2V1010) »W1010),V2010) ,U2(10) »W2(19) ,U3(10),V3(10), 


a j— weet et eZ =o ee ed 
© ww Re —— p~ wee >< 
~ < ee ele tet 2 be 
(om ~ — eat i Rape ates et 
Wc NS ee 3e00 © By eS Se | 
mc = He mnt om 
Ome met C) — 3. r AK 2 €or 
O~ oO f= tO. — PE a) + 0. am + See. et 
SS —t « « ~< «ke —~ + + 
Zz a2-oO o-~r HiN~ ~ = aT ~ =~ = 2 
~OW) W er « —~O it >oaxnx Cy CY —_—a~S O~ Im wt ee ee ee ee KE ome e 
— AQAAz- NQeewlU Om m<I | I Me —F- OL eb Oe COON 
s mee © OC -- LU eC bs Oe i wet Cet be re ee ee ed or Oe 1 ee 
al eS Oat to — Ore N SOS bh Celeron & Teer +e N 
\— m ODXSWNTO OxX~— ~< C195 wee OE ATOATOVE eG dE SUT Oe ee Cl LN 
UW LZ Ke TNO~O o& he et OX <{ ae (hy 2s cea ee) OO} tee O CE SoC wk tS at 
ea) mm Wee em om OTL O = wo eer KOO A) ae 8 A oe ete» OD aa | ipa! OR cee i tg aad 
— iz. Se wm we em OIL «+O QO »< ww yo wee we ei If STON et tO peng Ht bet et ee Yt 
az Samos ( ™ “Lh mm NS OR ew -~ -—J wa fod er ge Te mm we terete ee 
om eYOAaAIo erie S OWN SB eM MNM ISO MOE wit mK Mme bee KO OOK Keb rimx Le XX 
Mm o> OK TN Owe oe OE - a mtr (HE LLU eK ee ex ext rt oC perms et. oe <l meters Oe oe 
WW el NEKO KI LACK RO Oe OX = be ee, el CO~ Oz cee ee a Oe Se ere ta ~~ 
ZZ SF AS ree ee ew OIL Yee me ~— Denn) wirenKkK il UAL! wm R AAAs! COW ee K—OOAaAToORKO~ 
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\- ZMEZONNH~ OL OR 2 MON Ne NEN I mE NINN em HOO HI aati eeOD Fe KUTA ew KO Re oO 
5 NN SiN KE eo ON Se Oe we eT ESOC AM SCE KL KS 
DOM ™ eaqAQNOQANNOPOOCZ EADY ~O eer em~ew Er FEE OOOOH UWS KK NENINAMNS ~re NIH VLNGAN 
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APPENDIX B 


The load applied in the computer program is uniform and 
taken im the first mode only, i.e., n= 1. To accomplish 


this, the following changes were made to the program logic: 


SUBROUTINE ZANDZ 53 
DO 15 3 = 1, 4 659 
- IF(M.NE.1) ELLS(J) = 0 
RETURN 
END 
SUBROUTINE FORCE (Kk) 960 
DO 21 J3= 1, 4 081 
a IF(M.NE.1) ELLS(J) = 0 
RETURN 
END 
* = Required inserts 
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APPENDIX Cc 


The elements of the [A] and [B] matrices used in the 


eigenvalue solution are defined below. 


= z Zz ae ae 
Ay, = G1 7 F142 os) eee as) 2 
a 7 2 
Aig = Fi2% Rog = “E502 + Soo 
7 : 2 
#5 - *,3° oo 2S eee 
Fin = © Aog = Soa 
i 3 eo ia © 
M.. = G.. - E..0° A= G 
32 52 32 42 42 
a. -aem, = PwC ee Ne E..0° 
go 7 eso 33° 43 43 43 
_ F 2 - 
Az, = S34 7 £342 Aag = Saga 
B = 1/4 n? B = 1/4 nQ 
11 21 
: x 2 
B,, = 1/4 no Eee 17400 
a Se epee =—0 
Byy = 0 Bo, = 0 
B3, = 0 By, = 9 
a = 0 B,. = 0 
ee a 
Ba, °° Bie) ae 
Boy = © Pia = © 
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APPENDIX D 


(n) 


The following subroutine was used to calculate $. 
k 


and Gaee in the computer program. 


SUBROUTINE IMPERF 
Dinas lOn ENT (15) ,DFNT (15) 
COMMON /BLIMP/PHITB(200) ,PHIXB (200) 
AO=4 
SITGA=1000. 
EO=15200000. 
H=8./14. 
H1=1./H 
H2=1./ (2. *H) 
DOwl £1, 6 
AJ2=I+7 
BE EyN@weGO TO 8 
AJ2=0. 
S Peep (5,20) (PNP C1) gpI=8,15) 


f= 


O 
O 
NO 
A 
II 
NO 
bt 
n= 


2 DFNT (K)=(FNT (K+1) -FNT (K-1) ) *H2 
DIN? (1) = (PNT (2) -FNP (1) *H1 
DR (15)= (PNT (15) =FNT (14) ) *H1 
KR=I-1 
DO 3 K1l=1,15 
KZ=KR*15+K1 
PHITB (KZ) = (FNT (K1) *AJ2*EO/SIGA) * (0.00464) 
PHIXB (KZ) = (DFNT (K1) *EO/SIGA) * (-0.00464) 
3 CONTINUE 
i GONTINUE . 
MO FORMATV(GEI2.5) 
RETURN 
END 
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